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ABSTRACT 

A new method is presented to describe the evolution of the orbital-parameter dis- 
tributions for an initially universal binary population in star clusters by means of 
the currently largest existing library of iV-body models. It is demonstrated that a 
stellar-dynamical operator, £^y°,f ' rh (t), exists, which uniquely transforms an initial 
(t = 0) orbital parameter distribution function for binaries, T> ln , into a new distri- 
bution, Z) Mccl ' rh (i), depending on the initial cluster mass, M cc \, and half-mass ra- 
dius, r n , after some time t of dynamical evolution. For T> m the distribution func- 
tions derived by Kroupa (1995a, b) are used, which are consistent with constraints 
for pre-main sequence and Class I binary populations. Binaries with a lower energy 
and a higher reduced-mass are dissolved preferentially. The i7-operator can be used 
to efficiently calculate and predict binary properties in clusters and whole galaxies 
without the need for further iV-body computations. For the present set of TV-body 
models it is found that the binary populations change their properties on a crossing 
time-scale such that f^™ 1 ,rh (£) can be well parametrized as a function of the cluster 
density, p cc \. Furthermore it is shown that the binary- fraction in clusters with sim- 
ilar initial velocity dispersions follows the same evolutionary tracks as a function of 
the passed number of relaxation-times. Present-day observed binary populations in 
star clusters put constraints on their initial stellar densities, p cc \, which are found 
to be in the range 10 2 < p e d(^ r h)/-^o pc~ 3 < 2 x 10 5 for open clusters and a 
fewxlO 3 < p C ci(^ r h)/MQ pc -3 < 10 8 for globular clusters, respectively. 

Key words: star clusters: general - globular clusters: general - open clusters and 
associations: general - binaries: general - methods: N-body simulations - methods: 
analytical 



1 INTRODUCTION 

A significant fraction of stars in the sky are members of 
binaries or higher order multiple systems. The binary pro- 
portions thereby depend on the considered population. 

In the Galactic field (GF) about half of all late-type 
centre-of-mass (cm-)systems (where a system refers to either 
a sin gle star or a bina ry) are b in aries dp uqucnno y fc Mavoil 
19911: 



Mayo r et al 



Halbwachs et al.1 120031 ; 
20101 ). In contrast 



19921: iFischer fc Marcvl 1 19921: 
Raghavan et al. 20ld : iRastegaevI 



long-lived star clusters show a large 
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spread in the observed binary content. By analyzing the 
color distribution of main sequence stars, ISollima et al.l 
(2007) homogeneously estimated the global binary propor- 
tion for 13 low-density globular clusters (GCs, core-density 
log 10 p c /Mq pc~ 3 = —0.35 to 2.52) at high Galactic lati- 
tudes. Of all systems in GCs, between ~ 10 — 50 per cent are 
binaries. Their analysis revealed an anti-correlation between 
the binary proportion and clu ster age, that lies b etween 
6 and 12 Gyr in their sample. iMilone et al.l (|2008l ) had a 
larger sample and found an even stronger anti-correlation 
between the binary proportion and cluster luminosity 
(i.e. mass). A similar method was used to study five high 
latitude open cl usters with ages be tween w 0.3 — 4.3 Gyr 
homogeneously l|Sollima et al.l [2010l ). The binary fractions 
are generally larger than in GCs. In their cores a total 



2 Michael Marks, Pavel Kroupa and Seungkyung Oh 



of ~ 35 — 70 per cent of all systems are estimated to be 
binaries. Again, a dependence on cluster mass has been 
detected. 

Young star clusters (~a few Myr) show even higher 
binary fractions. The ~ 1 Myr old Orion Nebula Clus- 
ter (ONC) has a very high central d ensity (5 x 10 4 stars 
pc in the core, lMcCaughre"anll200ll ) and contains a pop- 
ulation whose binary proportion is comparable to that of 
the GF for the orbital period range P = 10 4 8 - 10 6 5 days 
llProsser et al.lll994l;|Petr et al.lll998l ) . A more recent analy- 
sis bv lReipurth et al.l (|2007n " shows a slight underabundance 
of visual binaries in the ONC in the semi-major axis range 
a — 67.5 — 675 AU. On the other hand, long period sys- 
tems (P = 10 7 — 10 8 ' 1 days) are signif icantly underrepre - 
sented in the ONC compared to the GF (jScallv et al.lll999h . 
The embedded cluster NGC 2024, also located in the Orion 
molecular cloud, has about the same age as the ONC but 
is less dense by a factor of about 10. The binary proportion 
there is significantly larger than in the ONC in the inter- 
val P = 10 5 - 7 - 10 7 1 days l|Levine et alj|2000l ). IC348 has a 
similar density and size as NGC 2024 but is 3-5 times older. 
The binary proportion in the range P = 10 5 — 10 7 9 days i s 
indistinguishable from that in the GF |Duchene et al.| [l999). 

Nearby low-density pre-main sequence T-Tauri star 
populations have a binary fraction up to twice as high as in 
the GF. Taurus- Auriga has a binary fraction of ~ 43 per cent 
in the separation range 18 - 1800 AU (1.9x the GF 
value) at an average stellar s urface density of a few stars 
pc~ 2 dKohler fc Leinertlll998h . The overall binary- fractio n 
in Taurus might be as high as 90 per cent ljDucheneHl999l ). 
The star forming region p Ophiuchus at about the same dis- 
tance as Taurus has a binary proportion of ~ 26 per cent 
in the separation r ange « 18 — 900 AU (1.1 X the GF value, 
iRatzka et aill2005l ). In the densest cores of p Ophiuchus cen- 
tral molecular cloud L1688 stellar den sities can be as high 
as 5000 stars pc~ 3 l)Allen et al] [20 02) . comparable to the 
ONC density, although in the ONC the high densities ex- 
tend over larger scales. The excess of binaries in Ch amaeleon 
comp ared to the GF is similar to that in Taurus l|Duchenel 
Il999h . The overabundance of binary systems is also seen for 
Class I protostellar objects (about 10 times younger than T- 
Taur i stars) in regions over the entire sky 1 Connellev et al.l 
l2008h . 

These findings could indicate that the star formation 
conditions in different populations have been different. How- 
ever, the binary proportion and properties in systems where 
interactions are important are generally not static. Consider 
therefore the relative equation of motion of a binary moving 
through a background of systems (e.g. in a star cluster), 
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Here, p — m\vnil(m,\ + 7712) is the reduced mass with mi 
and »ri2 being the masses of the primary and secondary com- 
ponent, respectively, e r is a unit vector pointing in the di- 
rection of the relative separation vector r, ui is the angular 
velocity of the binary and G is the gravitional constant. The 
relative motion may be perturbed by encounters with other 
systems, i.e. a time- varying additional acceleration, a pBrt (i), 
acts on the binary. Therefore binary properties such as the 
orbital period, P = P(i), the semi-major axis, a = a(t), 
and the mass-ratio, q = q(t) — 1112/111-1 < 1, become time- 



dependent due to stimulated evolution by interactions with 
cluster star 53 

The binding energy of a binary, being the sum of the in- 
stantaneous potential and kinetic energy of its components, 
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and thus depends only on mi, q and a, but not on the ec- 
centricity, e. Binaries can be classified according to Eb = 
— -Bbin > relative to the kinetic or thermal energy of an 
average cluster member, i?km = ma 2 c i/2, where m is the av- 
erage mass of a cluster star and a ec i is the three-dimensional 
velocity dispersion of the cluster. A hard binary is then re- 
ferred to as a system wich has Eb Skin and equivalently 
a soft binary has Eb -C -Ekin- 

An isolated self-bound stellar population generally seeks 
a state of energy equipartition, which leads to significant 
restructuring of the system when the two-body relaxation- 
time is shorter than its age. A hard binary therefore typically 
shares its energy during an encounter with the average clus- 
ter member, i.e. its internal relative velocity, v r , decreases 
and the components move towards each other. It follows 
that a and P decrease and thus Eb increases. A soft binary 
gains energy in an encounter (v r increases) such that a and 
P increase, and Eb decreases, eventually leading to dissolu- 
tion of the binary if the transferred energy is sufficient. This 
way, we qualitativ ely arrive at the Heqgie- H ills law of stim- 
ulated evolution (|Heggiel fl975l ; Ifflllsl il975l ): Hard binaries 
get harder and soft binaries get softer. Therefore hard bina- 
ries generate energy (energy sources heat the cluster) while 
soft binaries absorb energy (energy sinks cool the cluster, 
iKroupa. Petr fc McCaughreadll999h . 

It is for the above reasons that the binary proportion and 
properties can change if encounters are important, and an 
observed binary population can not be assumed to be primor- 
dial but may have instead been significantly larger at birth. 
In particular the GF distribution cannot be a ssumed to be 
the initial distributio n (Kroupa 1995a, see also lParker et all 
120091 ; iGoodwinllioiol ) . 

The evolution of binary popul ations in star clusters 
has b een studied throu g h iV-body dPortegies Zwart et al] 
1 19971 ; iTrenti et al.l 120071 ; IFregeau et al] I2009I. Oh et al. in 
prep.), Fokker-Planck and Monte-Carlo (iGao et al. 1991 



Giersz fc Spurzemll2000l ; [Fregeau et al. 2003 ; jyanoya et all 
2005) and analytical computations jSollimal 120081 ). high- 
lighting the importance of stimulated evolution. However, 
other processes alter a binary population as well. Bina- 
ries may form via ti dal capture, but this is generally 
an inefficient process jKroupal Il995al ; IBodenheimer et all 
1 19931 ; |l vanova et al. I 120051 ). Tidal-capture is possi ble to 
a larger extent in the dense cores of GCs l|Hut et al] 
Il992l ; iFregeau et al] 20091) only, but may al so often lead 
to mergers instead (IChernoff fc Huand fi996). Destruction 
of binaries may further occur via coalescence of compo- 
nents through t i dal dissipation between the components 
Hut et al. Il992l; Kr oupa 1995b) or hardening encounters 



Hilld ll984l ;~ lHurlev fc Sharal 120031 ). Unbinding of binaries 



1 The mass-ratio may change when a componen t of the binary is 
exchanged during strong encounters faillslll977h . 
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can happen through binary ste llar evolution, which is effi- 
cient in dissolving hard binaries |lvanova et alj|2005l ) . Direct 
collisi ons and mergers may occur in binary -binary interac- 
tions (jBacon et al.lll996l ; lFregeau et alj|2004 ). The contribu- 
tion of all these processes are however small in comparison 
with the destruction o f primordial or initial binaries through 
stimulated evolution (|Hut et al.lll992l ; ISollimall200SF ). 

The extraordinary high binary fractions in some young 
star clusters and for pre-main sequence T-Tauri stars are 
consistent with the assumption that all stars formed as bi- 
naries. Furth ermore, the vast majority of stars has formed in 
star c lusters (|Lada fc Ladall2003i ; lLadall2010l ; Brcs sert et al.l 
2010). In a series of three papers, iKroupal lfl995al lbllJl ex- 
plored the possibility that the different observed binary pro- 
portions might be the result of an environment-dependent 
dynamical evolution of an invariant, initially binary domi- 
nated population. In an attempt to understand the origin of 
the GF binary distribution, a typical birth aggregate and an 
initial period distribution for binaries, that lead after stim- 
ulated evolution to the observed GF distri bution when th e 
cluster has finally dissolved, was quantified (|Kroupalll995aT ). 
This Inverse Dynamical Population Synthesis shows that 
the GF population comes from short-lived, initially binary- 
dominated aggregates, which are dynamically-equivalent to 
the dominant mode cluster (Nb = 200 binaries and half-mass 
radius rj, = 0.8 pcjj In fi nding a solution to Inverse Dynam- 
ical Population S vnt hesis . I Kroupal (|l995bl ) derived the birth 
and the initial (or primordial) distribution for the periods, 
mass-ratios and eccentricities of binary stars. The method to 
find the birth distributions assumes that masses for stars are 
select ed from the stellar initial mass function (IMF. IKroupal 
120011 ). which is a two-part power-law in the stellar regime 
(£(m) oc m" Q \ ai = 1.3 for 0.08 sC m 0.5M o , a 2 = 2.3 
for m > O.5M0). The binary components are then paired 
randomly. Eccentricities are selected from a thermalized dis- 
tribution and periods are drawn from a distribution which 
first rises with increasing period and becomes flat for the 
longest periods. The binaries with the so selected properties 
then undergo a phase of re-distribution of energy- and angu- 
lar momentum within their circumstellar material, termed 
pre-main sequence eigenevolution designed to introduce cor- 
relations between orbital-parameters as observed, resulting 
in the initial distributions (Sec. 12. 2| ). 

In this paper a novel description of the evolution of bi- 
nary properties due to stimulated evolution is developed, 
which are init ially binary-dominated and start with the 
IKroupal l|l995bT ) initial binary distribution functions. We ap- 
ply this model to investigate the early binary evolution in 
A-body computations of star clusters by Oh et al. (in prep.) 
starting with different initial conditions (mass and size). The 
here performed analysis results in a tool to study, predict 
and compare binary populations and their orbital parameter 
distributions with observations without the need for further, 
computationally expensive, iV-body integrations. 



2 BINARY DISTRIBUTIONS 

The forthcoming sections deal with binary distribution func- 
tions (BDFs). In this section we will therefore introduce 
these. 

A population of binaries with primary-component 
mass mi (e.g. for G-dwarf binaries, mi m 1M©) is de- 
scribed by the distribution of their dynamical properties, 
2?(log 10 P,q,e : mi), i.e. the distribution of periods, P (mea- 
sured in days), mass-ratios, q, and eccentricities, e. For sim- 
plicity we assume that the orbital-parameter distribution 
can be separated, 



^(logio p > e > 9 : m i) = $io gl0 p(w,i) $,(mi) $ t (mi) 



(3) 



i.e. the quantities are independent of each other (which is 
true for P > 10 3 d binaries, Sec. 12. 2|) . We refer to <&i O g 10 p, 
$q and $ e as the period, mass-ratio and eccentricity BDF, 
respectively. Most commonly used among observers is the 
period BDF, since P is relatively easily accessible. 
The fraction of binaries with period log 10 P is 



Mio gl0 p) = 



N b Qog 10 P) 
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where iV;,(log 10 P) is the number of binaries with period 
logi P in an interval [log 10 P, \og 10 P + A log 10 P] and 
iVcms = A s , tot + Af,.tot is the total number of centre-of-mass 
(cms) systems in the population, i.e. the sum of all single 
stars and all binaries. The distribution of /& normalised to 
the width, A log 10 P, of a period bin, is defined to be the 
period BDF, 
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A cma Alog 10 P 
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Using this normalisation the area under the period BDF 
yields the total binary fraction of the population, 

Hog 10 P max 

/6,tot = / $i ogl0 p(mi) dlog 10 P . (6) 

" , log 10 P m i„ 

Here, log 10 P m i n and log 10 P max are, respectively, lower and 
upper bounds to the period BDF. As an example the period 
BDF of G-dwarf binaries in the GF (|Duquennov fc Mayor! 
119911 ) is depicted in the left panel of Fig. [1] as the solid 
histogram (Sec. I2.l"|l . The mass-ratio and eccentricity BDF 
(linear scale) can be defined in a similar way. 

For later analysis, and because it is assumed to be the 
more physical quantity, instead of $io Slo p(mi) this work 
will mainly make use of the energy BDF, $i ogl0 E b (mi). The 
birth and initial energy BDFs follow from the corresponding 
period BDFs (Sec. 12.21 using a Monte-Carlo method follow- 
ing the procedure in lKupper. Kroupa fc Baumgardtl (2008). 
For the A-body computations used here (Sec. [3} the energy 
for each individual binary in these integrations can be di- 
rectly calculated and the energy BDF is easily constructed. 



2.1 The Galactic field 



2 A cluster being dynamically equivalent to a different one is a 
cluster that evolves its binary population in the same way. 



In a long-term radial-velocity survey of nearby solar-type 
stars, the GF distribution o f binary properties has thor - 
oughly been investigated by IDuquennov fc Mavoij (|l99l[ ). 
Combined with data on visual binaries and common proper 
motion systems this survey found that the proportion of G- 
dwarf binaries (binaries which have a G-type primary) in the 




Figure 1. The left and right panels show the adopted period and energy BDFs, respectively. Both panels depict the same birth (weak 
solid lines, eq,|7|l. initial (thick solid lines, eq. [7] + pre-main sequence eigenevolution) and G-dwarf GF (solid histograms and dotted 
lines, iDuquennov h M ayor 1991) distributions. The energy BDFs in panel (b) follow fro m the period BDFs in panel (a) by applying a 
Monte-Carlo method by sampling binaries from analytical distribution functions (Sec. 12.21 K tipper et alj|2008h . For comparison with the 
birth and initial period BDFs, the symbols with errorbars in panel (a) show resu l ts from pre-main sequen ce observation of Taurus-Auriga 
(labelled can.PMS, f romlMathieulll994l; iLeinert et al1ll993l; iRichichi et al.lll994 iKohler fc Leinertlll998T). Lupus (Koehler , priv. comm.), 
Upper Sco A ('UScA. lBrandner fc Koehlerlll998l : lKohler et al.ll2000h and Class I protostellar objects l[Connellev et aLll200gt ). The GF BDF 
originates from the birth distribution, ©birth: after pre-main sequence eigenevolution and stimulated evolution in the dominant-mode 
cluster (Sec. [272J) . The eigenevolution operator, Cee (eq. (Tojl, transforms ©birth into the initial, i.e. eigenevolved birth distribution, X>; n 
fSec. I2.3")l . The GF distribution, T> Maol ' rh (t) = PqFi results from T> ln after applying the stellar dynamical operator, il^J^' rh (t) (eq. EJl, 
for the dominant-mode cluster (M oc i/Mq = 128, r^/pc = 0.8, t = 1 Gyr, Sec. 12 . 3 1 - 



GF is G /i,,tot « 0.53 ± 0.08. The distribution of orbital peri- 
ods (Fig. [T] solid histograms) is rather well approximated by 
a Gaussian distribution in log 10 P with mean log 10 P = 4.8 
and dispersion a\ i , n p = 2.3 (Fig. [T] thin dott ed curves). 
iMavor et all l| 19921 ) and lFischer fc Marcvl l| 19921 ) did a sim- 
ilar analysis for K- and M-dwarf binaries and found period 
distributions virtually indistinguishable from the G-dwarf 
distribution, i.e. the period BDF in the GF does not sig- 
nificantly depend on spectral type. The total proportion 
of late-type binary systems in the Galactic disc amounts 
to G ' K ' M f b tot — 0.47 ± 0.05, which is a weighted average 
of the proportion of binaries among, G, K and M dwarfs, 
respectively ( G / b ,tot = 0.53 ± ^ 0.08, K f b . to t = 0.45 ± 0.07, 
M /b,tot = 0.42 ± 0.09, lKroupa|[l995al l. 

The question of the origin of the GF binary population 
is debated. Since the binary proportion in various systems 
lies between 10 < /f,,tot < 90 per cent (Sec.[T]) and is in many 
cases different from the GF population, authors hav e argued 
for e nvironment-dependent binary formation (e.g. iKroupal 
l201ll . and references therein). 

The picture of binary formation may however be unified 
if one assumes that all stars form in small clusters and as 
members of binaries or higher-level multiple systems which 
are subsequently partially removed due to energy transfer in 
encounters between sy stems in stellar populations, which is 
the ansatz followed by IKroupal |l995al lbllcT) . This is the topic 
of the next subsection. 



2.2 The birth & initial BDF 

If all star formation takes place in small-sized aggregates 
with a few hundred late-typ e stars, all of them initially con- 
tained in a binary, IKroupal (|l995al ) noticed that it is pos- 
sible to find a typical birth configuration, that leads to the 
same BDFs as that of the GF after aggregate dissolution. 
He finds that a dominant mode cluster exists from which 
the typical GF binary originates. This cluster has initially 
fVcms = A'&.tot = 200 (i.e. f b>tot = 1) binaries and a half- 
mass radius of rh = 0.8 pc and evolves its binary po pulation 
to resem ble that of the GF. Inverting the problem. I Kroupal 
(|l995al ~) was able to infer a period BDF which turns into the 
GF period BDF after dynamical evolution in the dominant 
mode cluster (Inverse Dynamical Population Synthesis). 

Selecting binaries from the so derived period BDF, ran- 
domly assigning component masses from the stellar IMF and 
selecting eccen tricities from a thermalized distribution (e.g. 
IKroupal l2008bl . /b(e) = 2e), does, however, not lead to the 
observed correlated distributions between mass-ratios, pe- 
riods and eccentricit ies for main se quence systems in the 
GF with P < 10 3 d jMathieulll994ri Short- period binaries 
do not show eccentric orbits. IKroupal (|l995bl ) attributes the 
presence of such correlations to pre-main sequence eigenevo- 
lution (EE), which is the evolution of orbital parameters 
due to a re-distribution of energy and angular-momentum 
of short-period protobinary systems within their circumstel- 
lar discs and due to tidal dissipation and mass-growth dur- 
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Figure 2. The initial distribution of binding energies (top panels), mass-ratios (middle panels) and eccentricities (lower panels) for 
primaries of all masses before (left panels) and after EE (right panels) for different bins in period. The area under each distribution 
equals the total binary fraction in the corresponding period bin. The sum over the binned BDFs in each panel results in the corresponding 
full BDF. Top panels: Strongly bound, i.e. large-E;, binaries have a short period. Orbits with energies log 10 > 3 occur only after 
EE, when periods < 10 days occur (in all the left panels there is no continuous-line histogram for the shortest-period bin; compare also 
the thin with thick solid line in Fig. [T^). Middle panels: The component masses of all birth binaries are selected randomly from the 
stellar IMF giving rise to the shape of the cj-BDF before EE. After EE short period binaries have a tendency toward high mass-ratios. 
In fact, all binaries with log 10 P < 1, which are not existent before EE, have a mass-ratio close to unity. A peak at high-cj often seen 
in observations thus occurs. Lower panels: The eccentricities of all birth binaries are selected from a thermalized distribution (left). EE 
tends to circularize orbits (e — > 0) which had a high eccentricity at birth (since these have small peri-astron distances). All the newly 
formed log 10 P < 1 orbits are essentially circularized. After EE, orbits with log 10 P < 3 have a bell-shaped e-BDF and are increasing 
otherwise, as in the G-dwarf data of DM91. 



ing the stellar proto- and early pre-main sequence stellar 
evolution phase. The uncorrelated q-, P- and e-distributions 
of a protobinary birth BDF with randomly sampled stellar 
masses from the IMF and which are subject to pre-main se- 
quence EE turn into a new initial or eigenevolved birth BDF 
which shows the observed correlations. His model suggests 
a birth period BDF for late-type protobinary systems of the 
form 



C G p("H<2M0) = 2.5 



logic P ~ 1 
45 + (log 10 P - if 



(T) 



which has a maximum period log 10 P max = 8.43 (Fig. H]). 
Binaries with periods P < 10 d do not occur in eq. [7] (thin 



solid line in Fig. [T]). The model in iKroupal (|l995bl ) for pre- 
main sequence EE then evolves this birth pre-main sequence 
period BDF into the initial period BDF which has hard bi- 
naries ranging down to periods P ~ 1CP 1 d (thick solid 
line in Fig. consistent wi th pre-main sequence and Clas s 
I protostellar binary data. iKroupa fc Petr-Gotzensl l201ll ). 
Pre-main sequence EE acts in the same way in all star 
clusters and late-type stellar populations. Hereafter we will 
therefore refer to the birth and initial BDFs as the bi- 
nary orbital-parameter distribution functions before and af- 
ter pre-main sequence EE, respectively. The effect of EE is 
summarised in Fig. [2] 

Stimulated evolution of the initial, i.e. eigenevolved 
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birth population in the dominant mode cluster until aggre- 
gate dissolution after ~ 1 Gyr results in a population which 
has a comparable binary proportion and correlated distri- 
butions as the GF. 



2.2.1 Invariance of the birth & initial BDF 

The initial period BDF has not been derived by explicitly 
matching to observational pre-main sequence data, but by 
Inverse Dynamical Population Synthesis. Nevertheless, the 
outcome (eq. [7} matches the pre-main sequence and Class I 
protobinary data excellently (Fig.[]J,). The question may be 
raised if the initial period B DF is universal or if it varies be - 
tween star- forming systems (|Kroupa fc Petr-Gotzensll201lf ). 
Current star formation theory and hydrodynamical simula- 
tions have not been able to reproduce a distribution that 
resembles eq. (0. It is noted t hat the compu tations of iso- 
lated binary star formation by iFisherl (|2004T ) resulted in a 
log-normal period BDF broadly similar to that observed in 
the GF but he is not able to make specific predictions of the 
form of t he birth b inary f raction. 

The iKroupal l|l995bh birth and initial period BDFs 
agree with constraints from pre-main sequence data which 
have been obtained for distributed star formation and are 
rising towards long periods (Fig.lJj,). This suggests that the 
initial period BDF may be representative of star formation 
in general and may not be a strong function of the star 
formation envi ronment, although theo retical work suggests 
otherwise (e.g. iDurisen fc Sterzik|[l993 ). An invariant initial 
BDF would be quite similar to the unversality hypothesis 
for the stellar IMF which also follows from empirical evi- 
dence rather than theoretical considerations. In fact, since 
both distribution functions are the result of the same star 
formation process, and since the IMF is a result of processes 
"one-level deeper down" than the initial BDF, one can en- 
tertain the n otion that the initial period BDF ought also to 
be universal |Kroupall201l| y 

Note that /b,tot = 1 at birth independent of cluster 
density is a formal mathematical statement of invariance. 
Although wide binaries would not form in a dense clus- 
ter, they might originate from an initially much more ex- 
tended star forming region where /{,, t ot = 1 and the ini- 
tial BDF can be assumed before the cluster evolves into a 
denser configuration, which would be the case when the clus- 
ter forms dynamically cold (| Walsh et al . 2004; Pc retto et al.l 
l200rj ; lLada et alj|2008t ) . 

In the foll owing, b i naries are assumed to be initially 
drawn from the Kroupa (1995b) initial period BDF (eq. [3 + 
pre-main sequence EE). Although it has, strictly speaking, 
been derived for late-type binaries, here the same distribu- 
tion for earlier types is used since no good constraints are 
available and since there is no currently significant indica- 
tion that it should be different for more massive binaries. 
Binaries with a massive primary (mi > 5Mq) pair up, how- 
ever, differently from low-mass binaries (Sec. [3j . 

2.3 Evolution of binary orbital-parameter 
distributions in star clusters 

The final goal of this investigation is to find an effective de- 
scription for the evolution of BDFs in initially binary domi- 



nated systems. The model to analyse the iV-body computa- 
tions below (Sec. |3j is devised here. 

Assume that a single star formation event produces 
an embedded clustei[f| with total mass, M ec i, in stars with 
Nb(t — 0) binaries in a region which has half-mass ra- 
dius rh- The Nt(t = 0) binaries with primary-star mass 
mi have a distribution of initial dynamical properties, 
©i n (log 10 Eb, e, q : mi)Q. After some time t of stimulated 
evolution in a cluster that is characterised by (M ec i, J"h), the 
distribution function of dynamical properties has changed 
to © M ° cl ' rh (t)(log 10 E b ,e,q : mi). The superscript (M ec i,rh) 
denotes that the resultant distribution depends on the start- 
ing, i.e. star-formation conditions while the initial stellar dis- 
tribution function might be invariant (Sec. 12 . 2 . 1 p . Since the 
mapping 

© in — ►©^'^(t) (8) 

is unique for each cluster (M ec i, ru), the st ellar-dynamica l 
operat or, f2 dy ' n =1 ' rh (t), can be introduced (|Kroupal 2002, 
1200881 ) such that 

pMeci.n, (t) = Q^Wh (i) a D in , ( 9 ) 

where g) is an operation(3 Given T>i n , eq. (J5J) can be 
solved by performing A r -body calculations. One solution 
for n dj ^ 1 ' rh (*) ^ s derived in Sec. [3] Note that for stim- 
ulated evolution in the dominant mode cluster 0,^y n = 

n 128M ,O.8pc (1 Qyr) _ 128M O. 8P C 1 ^ 

dynamical operator discovered by IKroupal (| 1995al fbh - 

Similar to the stellar dynamical operator one can also 
introduce an EE operator, Qee, that evolves the protobinary 
birth BDF, ©birth, into the initial, i.e. eigenevolved birth 
distribution: 

©in = S2 EE ® ©birth • (10) 

A solu tion for Qee has already been devised in IKroupal 
l|l995bl . see also Sec. l23)l . 

The full equation that transforms the proto-binary birth 
distribution into the dynamically evolved distribution then 
reads 

V M ^(t) = n^ r » (t) ® ^EE ® ©birth , (11) 

which captures both the effects of pre-main sequence EE 
and stimulated evolution. 

3 iV-BODY INTEGRATIONS 

In order to find a solution for ^yn' '' h (*)> we utilize N- 
body integrations of star clusters performed by Oh et al. 
(in prep.). Therefore the main features of these integrations 

3 An "embedd ed cluster" needs not necessarily to evolve to a 
bound cluster llBoilv fc Kroupa 2003a]|3) and is taken to be any 
group of freshly formed stars with a surface density >few stars 
pc -2 . 

4 Note that we use the binding energy instead of period since we 
will perform our analysis of the models in energy space. 

5 For example, considering for simplicity only the initial period 
BDF, it can be described by an array of 3>i o g 10 p values, the op- 
eration is then a matrix multiplication yielding the new $i ogl0 p 

values. The matrix, fijJ5f ,ri * (t), is diagonal. 
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Table 1. Initial and final properties of the 18 computed models. 
Each integration of a type [M = M ec i/M@,rh/pc) was repeated 
100 times with a different initial random number seed and the 
results presented in Sec. I3,2l are averages over all such renditions. 
The following four columns give the initial values for the mass- 
density within the half-mass radius (in Mq pc -3 , eq. I13H . the 
crossing-time (in Myr, eg. 1121) . the median two-body relaxation 
time (in Myr, eg. I14|l and the velocity dispersion (in km 2 s -2 , 
eg. 1 15D ■ The last two columns show the half- mass radius (in pc) 
and the binary fraction after 5 Myr of evolution. 







initial values 






final 




lg(M) 


Hp) 




*rel 


ecl 




fb 


0.1 


1 


3.08 


0.30 


0.23 


0.17 


0.40 


0.54 


0.1 


1.5 


3.58 


0.17 


0.27 


0.55 


0.76 


0.50 


0.1 


2 


4.08 


0.09 


0.36 


1.74 


0.92 


0.44 


0.1 


2.5 


4.58 


0.05 


0.51 


5.51 


0.67 


0.40 


0.1 


3 


5.08 


0.03 


0.76 


17.42 


0.65 


0.35 


0.1 


3.5 


5.58 


0.02 


1.16 


55.10 


0.60 


0.29 


0.3 


1 


1.65 


1.56 


1.18 


0.06 


0.30 


0.80 


0.3 


1.5 


2.15 


0.88 


1.40 


0.18 


0.36 


0.74 


0.3 


2 


2.65 


0.49 


1.87 


0.58 


0.54 


0.65 


0.3 


2.5 


3.15 


0.28 


2.66 


1.84 


0.64 


0.57 


0.3 


3 


3.65 


0.16 


3.95 


5.81 


0.94 


0.52 


0.3 


3.5 


4.15 


0.09 


6.02 


18.37 


0.90 


0.45 


0.8 


1 


0.37 


6.79 


5.16 


0.02 


0.64 


0.95 


0.8 


1.5 


0.87 


3.82 


6.12 


0.07 


0.75 


0.94 


0.8 


2 


1.37 


2.15 


8.16 


0.22 


0.75 


0.90 


0.8 


2.5 


1.87 


1.21 


11.60 


0.69 


0.82 


0.82 


0.8 


3 


2.37 


0.68 


17.20 


2.18 


1.07 


0.74 


0.8 


3.5 


2.87 


0.38 


26.21 


6.89 


1.13 


0.65 



are re-called and the results necessary for our analysis are 
presented before we proceed to the analytical description. 



3.1 Setup 

Each integration starts with a given total mass, M cc \, in 
stars and a half-mass radius, r h . The grid of models com- 
puted together with some basic quantities is summarized 
in Tab. [T] The computations cover a range of initial den- 
sities, 5 < p cc i(^ r h ) < 7.55 x 1O 5 M pc -3 . An exten- 
sion to more massive star- forming events, M ec i > 1O 3,5 M0, 
is currently computationally not feasible due to demand 
on CPU time by the binary-rich clusters used here. Stel- 
lar masses for N* st ars in the integ rations are selected from 
the canonical IMF l|Kroupall200]] ) between the hydrogen- 
burning mass-limit of 0.08M Q and the maximum stellar 
mass m max deriv ed from the maximum ste llar mass - clus- 
ter mass relation (jWeidner fc Kroupall2006J ). The mean mass 
thus corresponds to the average mass of the stellar IMF, 
fn ~ 0.4Mq. All stars are members of a binary initially 
[fb,tot(t — 0) = 1] whose period is selected from the birth 
period BDF (eq. [7|. Eccentricities are drawn from a ther- 
malize d distribution w hich is not affected by stimulated evo- 
lution l|Kroupalll995bh . Stars with masses m* < 5M Q are 
paired randomly to form a binary (random pairing), while 
stars more massive than 5Mq are first sorted by decreas- 
ing mass and then paired with the next massive one so that 
massive binaries likely have high mass-ratios (ordered pair- 
ing). This procedure is chosen to mimick reality, since ob- 
servatio ns indicate that massive stars prefer a massive com - 
panion (|Kobulnickv fc Frverl [20071 ; ISana et al.|[20ol . 120091 ). 
Masses, periods and eccentricities are subsequently changed 



according to the pre-main sequence EE recipe (Sec. 12.2] ) be- 
fore the simulations are started in order to have realistic ini- 
tial conditions. Additionaly, the code adjusts the semi-major 
axes of binaries with the closest peri-astron distances before 
the integrations start to prevent physical contact between 
the components. This implemented procedure enlarges or- 
bits with periods of roughly log 10 P < until such con- 
tact doesn't occur any more. The centre-of-masses are dis- 
tributed according to a Plummer density profile initially and 
their velocities are chosen according to the correspondin g 



their velocities are chosen according to the corresponding 
velocity distribution (|Aarseth et alj 1l974l ; iKroupal BboSbf) 



The model clusters orbit a point-mass galactic potential at 
D = 8.5 kpc distance from the centre within which a total 
mass of Mg a i = 5 x 1O 1O M0 is contained. 

Due to the large number of models (100 integrations 
x 18 combinations = 1800 integrations in total), the com- 
putation time is restricted to the first 5 Myr of evolution. 
Therefore stellar evolution affects only the few most mas- 
sive stars exceeding « 4OM0 , which experience a significant 
change in their mass. The number fraction of stars more 
massive than 40Mq in the models is less than 0.1 per cent 
for a 1O 3 ' 5 M0 model, thus massive binaries hardly affect 
the overall binary population of the cluster. This on first 
sight short time-scale corresponds to 2> 10 dynamical (cross- 
ing) times for the densest configurations and it is argued 
that most of the binary evolution has already finished by 
then or proceeds only slowly, as oth e rs have done before 
jKroupal Il995al; iDuchene et alj 1 19991 ; iFregeau et ail 120091 : 
iParker et alj|2009h . 

Triple or hi gher multiplicity systems do not exist in the 
iKroupal fl995b) model. Higher-order multiple systems can- 
not be a significant contribution to the pre-main sequence 
population because if they were then their decay on a few 
crossing-times (<few 10 5 yr) would lead to a large frac- 
tion of single stars in star-forming regions which is not ob- 
served (|Goodwin fc Kroupal [2005). Higher-order multiples 
might instead form dynamically during the JV-body inte- 
grations and would be treated by the direct iV-body code. 
Such events are however unlikely and a higher-order mul- 
tiple system would likely not be long-lived (Sec.[T]). In the 
integrations, the higher-order multiplicity-fraction, 



Nt + Na 



/t + q 



N s + N b + N t + N q 



where Nt and N q are the numbers of triple and quadruple 
systems, respectively, does typically not exceed 1 per cent 
in individual time-steps. Therefore higher-order systems are 
not considered. Possible coalescences of binaries are treated 
by the code. 

Each particular (M ec i, rh) model is integrated 100 times 
with a different initial random number seed. All quanti- 
ties stated in the upcoming sections are averages over all 
such computations. This procedure is well suited for the 
present purpose since the aim is to investigate the typical be- 
haviour of a cluster with a given initial condition (M cc i, r^). 
The influence of the tidal-field is small, i.e. the vast ma- 
jority of stars stay within the tidal-radius during the 5 Myr 
integration-time (Sec. l3.2?2"j) . For the calculation of the quan- 
tities, at each time all systems in the integration are consid- 
ered and no systems are removed. 




Figure 3. Time-evolution of the binary-fraction in some model clusters in dependence of (a) the number of crossing-times and (b) the 
number of relaxation-times. Models having the same density (t cr oc p~^' 5 ) in panel (a) evolve their binary- fraction in the same way, while 
the corresponding tracks in panel (b) lie apart. Models with similar initial velocity dispersion but different i cr and t re i follow similar 
tracks in panel (b). In panel (a) the numbers are log 10 (p e cl/^0 P c_ 3 ) ana m panel (b) the numbers are ff^ d /km 2 s~ 2 . 



3.2 Results 

3.2.1 Time-scale for binary evolution 

The question to be answered here is on which time-scale a 
population of binaries changes its properties. In stellar dy- 
namics two ti me-scales are important. T he first one is the 
crossing-time ijBinnev fc Tremaind 120081 . using the virial- 
theorem and appropriate units), 



TT~ = 3-0 ( 

Myr 

where 



100M 



V M, 



ccl 



Peel = 



3M CC 1 



(12) 



(13) 



stands for the mass-density within the half-mass radius. The 
second important tim e-scale is the half-mass or median two- 
body relaxation-time l|Spitzer fc Hartlll97lh . 



21 



1M / M ccl 



Myr ln(0.4 X N) m \ 100M© 



(14) 



where N is the number of stars in the population. The 
crossing-time is the time which it takes a typical cluster 
member (v* ~ a, m* ~ m) to orbit the cluster once at 
a given radius, while a significant re-distribution of energy 
occurs on a two-body relaxation time-scale. 

Fig. [3|a) shows that in clusters with approximately the 
same density (i.e. the same t C r, eg. I12[) the binary-fraction 
decreases in a self-similar way besides having more than 
one order of magnitude difference in mass and more than 
a factor of 10 difference in their respective relaxation-times 
(Tab. [I}. Generally, the larger the stellar mass-density, the 
lower is the resultant binary-fraction. It can be seen, how- 
ever, that among the models with comparable density the 



one with the slightly lower density shows unexpectedly also 
the lower binary-fraction although this difference amounts 
to 2 — 3 per cent only (see also Fig. 1 101 below). Whether this 
is a result of, e.g., the averaging technique, due to stellar 
evolution (there are fewer massive-stars in lower-iV mod- 
els), or a real 2 nd -order effect is so far difficult to establish. 
An obvious difference between the models of similar density 
is that the ones with the slightly lower density contain a 
significantly lower number of stars, which might indicate a 
connection with relaxation processes (eq.QT}. However, since 
the difference occurs immediately after the integrations have 
started and afterwards their respective binary-fractions de- 
crease in the same way (the curves run parallel) it is likely 
that the first disruptions of systems occur quicker in the 
lower- N models simply due to the more compact configura- 
tion (smaller r^). 

In contrast, Fig. EJb) shows that the corresponding 
models with comparable densities evolve completely differ- 
ently on a relaxation time-scale. It can however be seen 
that, intriguingly, some models appear to follow comparable 
tracks in the sense that after a given number of relaxation- 
times they have the same binary-fraction. We note that the 
crossing-time as well as the relaxation-time in these models 
are completely different (Tab. Q]). It turns out that models 
which follow the same track have a similar initial velocity 
dispersion, 



2 

°"ccl 



>GM CC 
2r h 



(15) 



where s ~ 0.88 is a s tructure factor valid for a Plummer 
model l |Kroupj|20u8ah . The velocity dispersion in a clus- 
ter determines the current boundary between hard and soft 
binaries (Eb ~ -Ekin with .Ekin oc o~% cl , Sec. [TJ, such that ini- 
tially the same fraction of soft binaries is present which are 
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easily dissolved. This may be related to the velocity disper- 
sion in a cluster evolving on the energy-equipartition time- 
scale, i.e. clusters with similar velocity will, at the same dy- 
namical age (the same number of passed relaxation times), 
have broken approximately the same fraction of soft bina- 
ries. If true then this suggests a close coupling of the binary 
population and its host cluster. However, it is noted that the 
early evolution of the velocity dispersion may be driven by 
stellar evolution and binary-burning rather than relaxation 
effects. 

In order to understand the origin of this result one has 
to retreat to new iV-body integrations to disentangle effects 
driven by dynamics, stellar evolution and the presence of 
a mass-spectrum. Such computations are planned for the 
future, but exceed the scope of this analysis. Here, the aim 
is to devise an analytical description of the evolution of a 
binary population in a star cluster which implicitely includes 
all these effects. 

It is thus found that t cr is the time-scale over which the 
binary populations in the present low-mass clusters evolve 
through stimulated evolution, i.e. the cluster density is the 
primary parameter driving the depletion of the binary pop- 
ulation in these. 

A theoretical justification for this finding is however 
difficult. It is intuitive to assume that the denser a system is, 
the more encounters a binary experiences per unit time. But 
consider the number of collisions experienced by a binary on 
a circular orbit at the half-mass radius after one t cr . Let E co n 
be the cross-section for a single collision and pi oc the local 
density at rh, then, 



A^coii = 27rrh£ co npi oc m 



(16) 



For pi oc w p cc i, the number of collisions per unit time then 
becomes 



iVcoii / t CI oc r h p c ' c l 



3/2 



(17) 



In two populations of the same density but different sizes, 
a binary in the larger one will experience a larger num- 
ber of encounters since a longer path is travelled within 
the same time. Additionaly, clusters having different veloc- 
ity dispersions (eq. I15|) but the same number of collisions 
lead to different probabilities to break a binary. However, if 
a mass-radius relation for young clusters is absent or just 
very shallow dZepf et al|ll999l ; iLarsenl [20o3 ; iKroupal [20051 ; 
IScheepmaker et al.l l2007f ). this leaves p e ci as the parame- 
ter determining the number of collisions. But this argument 
does not hold for the models with diffe rent radii used here. 
Using an energy diffusion argument, iBinnev fc Tremainel 
(2008) show that the time-scale for binary destruction both 
through instantaneous ionisation through a single encounter 
and via multiple encounters, respectively, is oc Ccdp^d -1 . 
In any case, both strategies do not account for the rapid 
change of the half-mass radius and density over a crossing- 
time in the present models (Sec. I3.2.2p so that these ap- 
proaches might be too simple. 



3.2.2 Mass-, size- and density- evolution 

Since at all times all stars in the models are considered to 
extract the quantities of interest, none being removed, the 
mass of the population stays constant except for the negli- 



gible stellar evolutionary effe cts (Sec. 13.11) . Considering the 
mass within the tidal-radius l|Spitzerlll987h . 

only, the mass-loss for the most massive model clusters 
amounts to maximally 11 per cent of the initial mass and 
less than 1 per cent in the lowest mass clusters (including 
loss of systems over the tidal boundary + stellar evolution). 
The effect on the considered binary properties (Sec. I3.2.3[) 
is negligible: The difference in the binary-fraction consider- 
ing all systems or within the tidal-radius only is at most 
2 per cent for the present mass and size range. 

The size of the model clusters as measured by the 
half-mass radius rh (Fig. HI top left panel) stays generally 
constant for some time until the clusters start expanding. 
Standard secular evolution (no binaries) cannot be respon- 
sible for this expansion since the increase in half-mass ra- 
dius would be slower: A single-mass cluster with no pri- 
mordial binaries exh ibits a constant ha l f-mass radius over a 
much longer perio d dHeggie et al.ll2006l: iKupper et al.ll2008l : 
iHurlev fc Mackevll2010l : IConverse fc StahlerlbOllh . slow ex- 
pansion occuring only when binaries have formed dynami- 
cally, halt core-collapse and re-heat the cluster. In a multi- 
mass cluster without initial binaries a short period of radii- 
shrinkage is evident, befo re all mass-shells start expanding 
l|Converse fe Stahlerll201lT h The presence of binaries speeds- 
up the expansion of the half-mass radius considerably, as 
shown by the compute d isolated single-ma ss clusters with 
10 per cent binaries bv lHeggie et a l. (2006) (a factor 4 ear- 
lier than in the single-star, single- mass case). Adding a tidal- 
field, the removal of systems over the tidal-boundary com- 
pensates for the heating such that the half-mass radius stays 
constant and might even d ecrease, if the tidal-field is strong 
enough |Trenti et al.ll2007f ). 

Therefore the increase in size is due to the presence 
of a mass-spectrum which drives mass-segregation and the 
energy generation through hard binaries. The onset of ex- 
pansion occurs earlier and the expansion-rate is larger the 
denser the cluster is, since binary-burning is more efficient 
initially. Therefore compact configurations exceed more ex- 
tended configurations in size after some time of evolution, 
thereby slowing down binary disruption at later times in the 
initially densest clusters (Sec. [5}. 

The rate of expansion appears to be comparable for 
initially same-sized clusters when plotted as a function of the 
number of relaxation-times (Fig. [4] top right panel) while the 
onset of expansion happens at a lower number of passed-by 
t re i for higher-mass clusters. The initial t IC \ for higher- M BC \ 
models is larger than for lower-mass clusters of the same 
size. 

The second row in Fig.|4]shows the evolution of the core- 
to half-mass radius ra tio, r c /rh, r c being ident ified from the 
integrations following ICasertano fc Hutl (| 19851 ). The cluster 
core starts to contract immediately while rh is still constant. 
This increases the stellar density and therefore the efficiency 
of binary disruption in the core. At later times the decrease 
of r c /rh is supported by an increasing rh. When r c /rh is de- 
picted as a function of time in units of t Te \ the models rather 
group to follow similar tracks as a function of their initial 
mass. Models with masses below 1OOA/0 should be looked 
at with caution since the determination of radii in individ- 
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Figure 4. Evolution of the half-mass radius (top panels), the core- to half-mass radius ratio (middle panels) and the density within 
the half-mass radius (lower panels) as a function of time in Myr (column a) and in units of the initial relaxation-time (column c). The 
line- types differentiate cluster masses, the line-strengths denote the initial half-mass radius. The denser a model is initially, the stronger 
is the expansion. Initially compact clusters can quickly exceed initially more extended models in size due to energy generated from hard 
binaries. Core-contraction is visible throughout the whole integration time, enhancing the stellar density in the core leading to more 
efficient binary disruption. The mass-density within the half-mass radius stays roughly constant for all initially spatially extended models 
(r^ = 0.8 pc) and decreases up to two orders of magnitude for the densest configurations within only 5 Myr. 
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ual realisations of a low-mass cluster is difficult given the 
low number of systems in these clusters (only 28(76) stars, 
or initially 14(38) binaries, in a model with 10(10 ls )Mo, 
respectively). 

In single-mass clusters with no binaries the initial core- 
collapse is rather de ep before gravothermal oscillations set in 
(Hcgg ie~et al.ll2006h . Introducing binaries, the core shrinks 
more rapidly (due to mass-segregation), but it is not so deep. 
After a few i re i, the ratio levels off and decreases only slowly 
over many relaxation-times. Thi s is in quali t ative agree- 
ment with the runs of r c /r h in iTrenti et al.1 (|201Ch . who 
computed multi-mass models with initial binary-fractions 
up to 10 per cent and without binaries. All their compu- 
tations show an initial decrease of the radius-ratio until 
they level-off to a more or less common value of r c /rh af- 
ter a few relaxation - times which i s of th e same order as in 
iHeggie et ail l|2006l ). Iliurlev et al.1 \200H ) computed models 
with 10 5 systems, a mass-spectrum and up to 10 per cent bi- 
naries, showing that after a phase of stellar evolution driven 
expansion r c /rh decreases continously over the computing 
time of 15 Gyr. Over the limited time-span in the present 
iV-body computations the ratio continues to decrease at the 
same rate, although the models reach « 10£ rc i ■ But this is 
rather driven by the continously rising half-mass radius than 
a contracti ng core, and r c /rh h as n ot yet reached the v alue 
reported in ITrenti et ail (|2010h and lHeggie et al.l l|2006l ). 

The lower panels of Fig. [3] show that the density within 
the half-mass radius, p oc i, stays more or less constant for the 
0.8 pc-sized models and decreases only slightly for M cc i < 
10 2 Mq and rh = 0.3 pc. For models in excess of 10 3 Mq with 
Th = 0.3 pc as well as for the low-mass rh = 0.1 pc model, 
p e ci decreases by more than one order of magnitude within 
only 5 Myr of evolution. For initially even denser models p oc i 
decreases by more than two orders of magnitude, thereby 
shrinking below the actual p cc i of initially much less dense 
clusters. No particular trend with t rE \ is seen in the density- 
evolution. 

3.2.3 Evolution of orbital-parameter BDFs 

The A^-body integrations show the dependence of the evo- 
lution of the binary population on the initial cluster param- 
eters (M G ci,?"h). In this section we describe the evolution of 
the energy, period, mass-ratio and eccentricity BDF in single 
clusters and between different models. The energy BDF will 
then be used to describe the binary evolution analytically. 
A more extended investigation of similar TV-body models 
and covering other aspects will be available in Oh et al. (in 
prep.). 

A population of binaries subject to stimulated evolution 
disrupts systems with the largest cross-section first. These 
are the binaries with the largest semi-major axes, a, longest 
periods, P, and, equivalently, the lowest energies, Ef,. The 
hard binary population is more or less unaffected by stimu- 
lated evolution (The Heggie-Hills law, Sec.[l|. This effect can 
be seen in the panels of Fig. [5] The models depicted in the 
panels of column (a) show the evolution in a cluster which 
starts with (10 2 ' 5 Mq,0.3 pc). Stimulated evolution causes 
the total binary fraction (the area under each distribution) 
to decrease. Low-_E(, (long-period) binaries are efficiently re- 
moved from the population and the peak of the distribution 
shifts to higher energies (shorter periods). After 1 Myr the 



binary fraction in this particular cluster shrinks to 71 per 
cent. At the end of the integrations (t — 5 Myr) the binary 
fraction is 57 per cent. The evolution of the binary propor- 
tion is therefore fastest at the start of the integration and 
the rate decreases at later times. The shape of the mass-ratio 
and eccentricity BDF is conserved with time. 

The models in the panels of column (b) of Fig. [S] 
compare the binary populations in three clusters with the 
same mass (M cc i = 1O 2 I5 M ), but different half- mass radii 
(rh = 0.1,0.3,0.8 pc) after t — 3 Myr of stimulated evolu- 
tion. It can be seen that stimulated evolution is most effec- 
tive in the densest configuration (rh = 0.1 pc). The denser 
the cluster initially the lower the binary fraction after the 
same time of evolution. 

The same trend is seen in the rh = 0.3 pc models with 
different stellar masses at t = 3 Myr (panels of column (c) 
in Fig. [5}. The more mass is inside the same volume, i.e. the 
denser the configuration is, the lower is the binary-fraction 
after the same time of evolution. 

We note that stimulated evolution can not account 
for the very low-Eb (long-P) binaries observed in the 
GF (Fig. [T]), since disr uption dominates over widening 
|Kroupa fc BurkertllioOlh . Thus, the low-energy part of the 
energy BDF is a forbidden region for binaries in a star 
cluster. While such binaries might exist, the life-time of 
these should be extremely short and the possibility of ob- 
serving such a system is low. A possible channel by which 
wide binaries with P > 10 9 d can exist in the field is 
their dynamical formation through co-moving stars that 
pair up after residual-gas e xpulsion from embedded clusters 
i|Kouwenhoven et alJl2oToh . 



4 ANALYTICAL DESCRIPTION 

This section seeks to find a solution for the stellar dynami- 
cal operator, fi^™ 1 ' 7 " 11 (t), by applying eq. (JSJ to the iV-body 
models. In the following we will distinguish between the N- 
body data, best-fit parameters/curves and fitting functions. 
Best-fit refers to fits to the set of individual iV-body in- 
tegrations with the same initial conditions (M ec i,rh) while 
the fitting-functions describe the behaviour of the best-fit 
parameters for all computations. The fitting-functions are 
used to calculate the model BDFs. 

4.1 Energy BDF 

Let r2dyn 1,rh (*) act on the ener gy BDF (via eq. [9jl, 

<£oE b (*) = <n"'' h (*) ® *log 10 E b ,in , (19) 

and the operation g) be a multiplication. Then fi dy ^ 1,rh (t) 
is simply the ratio of the resulting and initial energy BDF 
at some snapshot of the simulation. 

That way an empirical solution for fi t i™ 1 ' rh (t) can 
be obtained. The stellar-dynamical operator is a relatively 
steeply rising function with increasing binding energy which 
flattens at some level. Fig. [6] depicts the situation for a clus- 
ter with Mcci = 1O 3 M0 and rh = 0.3 pc initially. At any 
given time the position, height and steepness varies with 
the initial conditions: as the energy BDF of a less dense 
initial configuration peaks at lower energies and retains a 
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Figure 5. The evolution of the energy (top panels), period (second row), mass-ratio (third row) and eccentricity (bottom panels) BDF 
as function of time, t (column a), half-mass radius, r h (column b), and embedded cluster stellar mass, M ec i (column c). The fixed and 
varying parameters for each column are indicated in the top panels. Starting from the initial BDFs (solid lines and histograms, cq. [jj 
+ EE), stimulated evolution dissolves binaries and alters the BDFs. Wide binaries (low-i?;,, long-P) have the largest cross-section (a 
large semi-major axis) and are therefore removed first. The hard binary population more or less resembles the initial distribution (The 
Heggie-Hills law, Sec. [TJ. (a) As time progresses an increasing number of binaries from the low-energy (long-period) end are removed 
such that the peak of the distribution shifts to higher energies (shorter periods). The shapes of the mass-ratio and eccentricity BDF arc 
roughly preserved. Only the change in the binary fraction, i.e. in the area below the distributions, with time is visible. (b,c) The smaller 
or more massive, i.e. denser a cluster is, the more binaries are dissolved within the same time, as a result of higher encounter-rates in 
more compact configurations. 



larger binary population (Sec. I3.2[) . J7 d J£ 1 ' rh (t) also shifts 
to lower log 10 Et, the upper asymptote is closer to unity 
and the central part is somewhat steeper (see Fig. [7] below). 

The decrease of the empirically found operator (thin 
solid line) at the highest energies is a result of the initial en- 
ergy BDF used here which is found via a Monte-Carlo exper- 
iment (thick solid line, Sec. This trend is not seen in the 
pure A^-body based data (the thick grey line, arrived at by 
computing the ratio of the hist ograms in Fig. P i because the 
hardest binaries present in the lKroupal l|l995bl ) model have 
been moved to slightly lower energies before the integrations 
started (see Sec. I3.1[) . Furthermore the number of the hard- 



est binaries is not expected to decrease very strongly, only 
some may be ejected or d issolved in binary-binary interac- 
tions (|Fregeau et alj|2004 seen in the drop of the thick grey 
line) or moved to higher energy levels through hardening 
processes (seen by the grey solid line which exceeds unity in 
the last bin). The strong decrease is therefore not considered 
in the following. By using the Monte-Carlo based distribu- 
tion as the initial energy BDF (thick solid line) the operator 
never exceeds One as it always lies above the JV-body based 
distribution (thick grey histogram) for the hardest binaries. 
Some artificial hardening of binaries not present in the N- 
body integrations is however introduced due the neglection 
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Figure 6. The initial energy BDF (thick solid peaked line and the 
superposed thick histogram), <3?i O g 10 E b .im i s evolved into a new 



distribution, $ 



logl 



(i), after t = 3 Myr of evolution in a cluster 
with M ec i = 10 3 Mq and = 0.3 (thin histogram, left axis). The 
stellar-dynamical operator, f^y^ 1 '' 11 (t) (eq. 1191 thin solid line, 

right axis), that transforms 3>i O g 10 B 6 ,in — ► ^log *'^. W ' s f° un d 
by dividing the thin solid histogram (evolved BDF) by the thick 
solid line (initial Monte-Carlo based BDF) to be a rising func- 
tion which flattens at some level. The decrease at high energies 
is not considered since it is solely due to a difference between 
initial Monte-Carlo (thick solid line) and TV-body based energy 
BDF (thick histogram, Sec. 14. It . The pure TV-body based operator 
(thick increasing line) found by dividing the thin solid histogram 
(evolved BDF) by the thick solid histogram (initial TV-body BDF) 
doesn't show the decrease and hardening of the most strongly 
bound binaries is evident. The best-fit energy BDF (dashed lines) 
and the one resulting from the fitting-functions (eqs. I2U23I dot- 
ted lines) reproduce the TV-body computations reasonably well. 
Q^ ^ 1 ' 1 " 11 (t) acts as a transfer— or orbit- depletion— function (see 
also lKroupalll995ar ) . 



of the strong decreasing trend. The derived BDFs at the 
high-energy end are however hardly affected since at most a 
few hard binaries in a thousand are added that way. In par- 
ticular the stated points have no influence on the shape of 
the BDF at low-energies where the operator removes orbits 
efficiently and shapes the resulting BDFs. 

The empirical curves for fi d J£ !l ' Ph (i) can be described 
as the upper halves of S-shaped, or sigmoidal functions. A 
special case of a S-shaped curve is the logistic function, 



l(x) = 



1 + e- 



le[0: 1]. 



This sigmoid curve is centred around x = with l(x — 
0) = 1/2. For our purpose we take a generalized form of 
l(x) to represent the stellar-dynamical operator. The logistic 
function is rewritten to satisfy the requirements, i.e. 



dyn 



" h (£,t) 



A 



1 + exp [S (£ - £ cu t)} 



A 

2 ' 



(20) 



where £ = log 10 Eb is the x-coordinate of the energy BDF, 
A/2 = A MeaUrb (t) /2 < 1 is the level at which Q^'^ {£, t) 



flattens, S 

ness, or slope, and £ cut 



S eol ' rh (t) < is a measure of the steep 



£ 



(t) 



3glo(£b 



> 



—3.2 is identified with the logarithmic cut-off energy of 
the energy BDF below which no binaries exist, since Vt : 

<n" rh (£»t,t)=0. 

A, S and £ cut are free parameters that are determined 
for a fixed time and for each initial configuration separately. 
The procedure is as follows: Numerically, a grid over the 
relevant ranges of all three parameters is simultaneously 
covered and for each combination the corresponding stellar- 
dynamical operator (eq. I20j> is computed, which yields via 
eq. (|19[) the resulting energy BDF. The set of parameters 
(A, S , £ C ut) which minimizes the sum-squared difference be- 
tween the fitted and resulting TV-body energy BDF is then 



chosen to describe the wanted Q 



dyn 



(£,t). For the exam- 
ple in Fig. [6] the overlaid dashed lines represent the corre- 
sponding operator and energy BDF using the best-fitting 
parameters. The so determined set of values for A,S and 
£cut well reproduce the energy BDFs and binary-fractions 
of the integrations with differences of Aft — ±2 per cent in 
the worst case (Fig. 1101 below ). 

In order to be able to identify the correct stellar- 
dynamical operator from the initial conditions as identified 
by (Med, Hi) the dependence of A,S and £ cut on the den- 
sity is investigated. Fig. [TJa) shows that the cut-off energy, 
£cut, equals the low-energy boundary of the initial energy 
BDF for densities below w 1O 2 A7~0 pc -3 and then shifts 
to harder binaries with increasing time and density, reflect- 
ing that stimulated evolution depletes the binary population 
and is more effective in dense clusters (Sec. 13.21 eq. I12|l , 

The maximum, A, in eq. (|20[) decreases with increasing 
density (Fig. [7p) , again due to the more efficient depletion 
of binaries in high-density configurations. At earlier times 
the population is less evolved, so that ^4 is larger there. 

Fig. PT|f c) shows that the rising part is much steeper at 
lower densities, but is otherwise very similar for different 
times. This is understood since in low-density clusters £ cu t 
is close to the low-energy boundary of the initial energy 
BDF and in order not to remove too many binaries from the 
population, fi dy ™ 1,rh it) has to be steeply rising. S appears 
to flatten beyond 10 3 M Q pc -3 . 

The values of the parameters for different times become 
increasingly similar the higher the density is, showing that 
the binary-burning phase has already ended by 1 Myr for 
the densest configuration, since stimulated evolution is so 
efficient there. Additionaly it can be noted that the evolution 
slows down with time and should have more or less finished 
after 5 Myr for all densities (see also Sec. [3}. 

For all parameters we provide fitting-functions to de- 
scribe their behaviour in dependence of density and time. 

Fig. [7£a) suggests a linear increase with log p ec \ for the 
distribution of cut-off energies, where £ cut is different from 
the low-energy boundary. A linear decreasing function is 
suitable to describe the dependency on density of the height 
parameter, A, in Fig. [7jb) . Therefore the functions 



£cut 



a + &log 10 Peel 

-3.2 



if result > 
otherwise 



-3.2 



(21) 
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Figure 7. Values of the parameters describing the stellar-dynamical operator f2 dy "J 1 cl ' rh (t) (eg. 12011 in dependence of initial cluster density 
for different times. The fitting-functions are extrapolated to densities beyond the range of densities covered by the integrations which 
are denoted by the symbols, (a) The current (at time t) cut-off energy of the binary population, E cu t, increases with time and with 
increasing density, (b) The maximum, A, of the stellar-dynamical operator is lower at higher densities and later times, (c) The slope 
parameter, 5, takes larger values with increasing density (f2,^ 1 cl ' rh (t) gets flatter). The lines (solid: t = 1 Myr, dashed: t = 3 Myr, dotted: 
t = 5 Myr) depict the fitting-functions for the best- fit values of the iV-body data (Sec. 14. ip . All parameters become nearly identical for 
the highest- density computation. 



Table 2. Values of the coefficients in eqs. (J2TJ , (22j and {23]) 

that describe the dependency on the initial cluster density of the 
parameters £ cu t, A and S entering the stellar- dynamical operator 
(eq.HJ) for t = 1, 3 and 5 Myr. 





t 1 Myr 


1 


3 


5 




a 


-5.00 


-4.64 


-4.40 




b 


0.65 


0.59 


0.54 


A 


c 


2.15 


2.04 


1.97 




d 


-0.13 


-0.11 


-0.10 




e 


2.36 


1.24 


1.47 


s 


f 


2.34 


1.70 


1.35 




9 


0.89 


0.79 


0.82 



A = 



c + dlog 10 p cc i 
2 



if result ^ 2 
otherwise 



(22) 



where a, b, c, d are coefficients to be determined, are chosen 
to represent the data. 

The variation of the slope parameter with logp cc i 
(Fig. 0i) suggests a curve which is very steep at low densi- 
ties and flattens sharply at the highest densities. A function 
that fullfills these requirements has the form 



S = - 



1 



exp [e(log 10 pcci - /)] 



(23) 



where e, / and g are coefficients. 

The values for the coefficients of these fitting-functions 
for 1, 3 and 5 Myr of stimulated evolution are listed in Tab. [2] 
and the resulting functions are depicted in Fig.[7]as the solid, 
dashed and dotted lines. That these fitting- functions are able 
to reproduce the results of the A-body computations can 
be seen from a comparison of the binary-fractions with ft 
calculated making use of eqs. (|2Hl - (|23| l in Fig. 1101 below. 

Since the binary populations hardly change after 1 Myr 
of stimulated evolution in the densest configuration, we se- 



lect the parameters according to the t = 5 Myr fitting- 
functions beyond that density to extrapolate to higher den- 
sities, i.e. the parameters and therefore, the energy BDFs, 
are taken to be invariant for p cc i > 3.77 x 10 5 Mq pc~ 3 . 



4.2 P-, a-, e- and q-distributions 

A similar method is in principle applicable to calculate the 
evolved BDFs for P, a, e and q separately from their re- 
spective initial distributions. However, a different approach 
is used by noticing that the orbital-parameters are interre- 
lated by Keplers laws. 

Therefore at first a library consisting of ATiib binaries 
is compiled. Their periods are selected from the birth pe- 
riod distribution, their eccentricity is chosen from a ther- 
malized distribution. 2Nub component masses are randomly 
selected from the canonical IMF between the hydrogen- 
burning mass- limit (O.O8M0) and 150Mq. Stellar masses 
from this array are combined randomly if the primary mass 
is below 5Mq . Secondary masses for primaries more massive 
than 5Mq are selected from the array of 2Nm, masses such 
that birth mass-ratios are larger than 0.9 to mimick the or- 
dered pairing algorithm (Sec. [3}. EE is then applied to the 
so selected parameters to yield the initial binary properties, 
which are added to the library. The library is complemented 
by the mass-ratio, q, following from the eigenevolved com- 
ponent masses, the binding energy, 

'7rmini2 \ 2 / 3 



E h = 2 



-1/3 ^ 



and the semi-major axis, 
,mira2 



G- 



2E b 



(24) 



(25) 



In order to then extract the orbital-parameter BDFs 
using the known evolved energy BDF for a cluster of mass 
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Figure 8. Evolution of the P-, q- and e- distribution in the 
Af-body integrations (histograms) for the cluster with M ec i = 
10 3 Mq, rj, = 0.3 pc (same as in Fig. [H} and calculated as de- 
scribed in Sec. |4.2| (crosses) . The evolved distributions in the N- 
body integrations and the resulting BDFs from the analytical 
description agree very well. 



M ec i and half-mass radius (i.e. p ec i), £ c «t, A and S are cal- 
culated using the fitting-functions to determine Sl^^' 7 ^ (t) 
and the resulting energy BDF (eg. I19[) . From this immedi- 
ately follows the binary-fraction of the evolved population. 
The initial and evolved number of systems and binaries in a 
population is calculated from the number of stars in a clus- 
ter, JV» = M oc i/m (m w 0.4M Q for the canonical IMF), and 
from the binary fraction known from the initial and evolved 
energy BDF (eq. [6)l , respectively, 



N* 



1 + ft 



N b = /b.totiVc; 



(26) 



The iV; nit and N£ voWcd binaries are then distributed into 
bins according to the initial and evolved energy BDF, re- 
spectively, and the fraction of surviving binaries per bin, 



/surv — 



N mit 



(27) 



is determined. Of all binaries in the library, 1 — / sur v of them 
in the corresponding energy bin are replaced by two single 
stars. From the binaries left in the library the final BDFs 
for P, a, e and q are computed. 

Comparing the so constructed BDFs with the outcome 
of the JV-body computations for the Af cc i = 1O 3 M0, = 
0.3 pc cluster in Fig. [8] (same cluster as in Fig. [6]), very 
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Figure 9. Evolution of the BDF for the reduced-mass, fi (lower 
panel, histograms), of binaries in a small energy interval (—0.5 < 
log 10 Bf, < 0.5, see top panel) for the M ec \ = 10 3 ' 5 Mq, = 
0.1 pc model to test the ionisation probability on fi. The t = 
data represents the post-EE rt-BDF. The fraction of surviving 
binaries (eq. I27j l decreases with increasing fi (dashed line, right 
axis) as expected from eq. 1281 The initial number of binaries in 
the higher-/^ bins is, however, very low (see the text). The trend 
is qualitatively the same for all cluster masses but the initial fi- 
BDF for low-mass clusters (M ccl < IO^Mq) is confined to the 
first (few) bins since only low-mass stars are present. 



good agreement between the evolved model- and iV-body 
BDFs is found. This further suggests that the derived fitting- 
functions (eq. 1211 l22l and 123]) . dependent only on the initial 
mass-density, p ac \, are suitable to describe evolution of the 
binary properties in the present models since these formulae 
have implicitely been used to construct the evolved BDFs 
from the library. 



4.3 Reduced-mass dependent dissolution? 

The above procedure is correct only if the probability of 
destroying a binary is a function of binding energy onl y or 
other parameters have a negligible influence. iHutl (fl983T ) de- 
termined an ionisation cross-section, Ei on , through a large 
number of scattering events for the case of single— binary- 
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star scattering at high velocities which depends on the 
binding-energy and the reduced mass, n — mima/(mi + 
7712). His eq. (5.1') can be rewritten if the incoming single 
star (m.3, V3) is a typical cluster star (1713 as m, V3 as a ec i), 



Sion OC G 



2 /i m 



Eb a. 



eel 



(28) 



which suggests that the larger the reduced mass, the larger 
is also the probability of ionisation. Looking carfully into 
the iV-body calculations we indeed find Hut's ^-dependence 
when investigating the time-evolution of the /i-BDF (Fig. 0). 

The lower panel of Fig. [5] depicts the evolution of /i 
in a small energy interval where binaries are actively being 
burned (—0.5 < log 10 Eb < 0.5, upper panel). It is con- 
firmed that the larger /j, the lower is also the fraction of 
surviving binaries (dashed curve, right axis). Note that the 
trend is exactly opposite (more binaries survive in high-/i 
bins than in low-/i bins) when considering the full range of 
the energy-spectrum. The number of binaries in the highest- 
[s, bins is, however, much lower than in the low-/^ bins (< 3 
binaries per individual model cluster per bin above \i = 0.5 
for the highest-M ec i ) . The question is thus whether this 
/^-dependent dissolution has a significant influence on the 
BDFs constructed from the library, or not. 

Therefore one extreme is probed in which binaries in the 
library which are ionised and replaced by two single-stars 
are selected to be the ones with the largest [i. This method 
is expected to have the strongest influence on the BDFs 
compared to procedures in which binaries are destroyed ac- 
cording to the p-dependence in Fig. [§] In order to do so, 
the algorithm in Sec. 14.21 is followed, but all initial binaries 
within one energy bin are sorted according to their reduced- 
mass, such that those with the largest /x in a single bin are 
replaced first. The so constructed BDFs can be compared to 
the ones using the energy-criterion only. The resulting BDFs 
are indistinguishable from those in Fig. [5] (crosses), i.e no sig- 
nificant difference between a ^t-dependent and -independent 
construction method is found. 

Thus, although ionisation in the present calculations is 
dependent on the reduced-mass this has a negligible second 
order effect on the constructed BDFs. For this reason it is 
justified to remove binaries from the library as a function of 
their internal binding energy only. 



5 BINARY-FRACTION AND INITIAL 
CLUSTER DENSITIES 

Fig. [10] depicts the binary- fraction in dependence of the ini- 
tial cluster density for the iV-body computations, for the 
best-fitting parameters and for the analytical model (i.e. us- 
ing the fitting- functions) . The model description (Sec. [4} by 
means of the stellar-dynamical operator allows for an extrap- 
olation of the data beyond the maximum density reached by 
the iV-body calculations. The data suggest that for initial 
stellar densities as high as p cc i = 10 s Mq pc -3 the binary- 
fraction drops to 12 per cent. This has however not to be un- 
derstood as a lower limit to the binary-fraction in stellar sys- 
tems, but just to illustrate the overall disruption efficiency 
in star clusters of initially higher density clusters within the 
first 5 Myr. In fact, binary fractio ns of 5 — 7 per cent o nly 
have been reported for NGC 6397 (|Cool fc Boltonll2002l ). 
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Figure 10. Upper panel: Comparison of binary fractions in the 
TV-body integrations for three different times (solid lines, 1, 3 and 
5 Myr from top to bottom) and calculated for the best-fit pa- 
rameters (dashed lines) as well as for the model (i.e. making use 
of the fitting- functions, dotted lines). The model agrees very well 
with the computations. Beyond a density of safew xlO 5 M0 pc -3 
the extrapolated model predicts a further shrinking of the binary- 
fraction at even higher densities, reaching 12 per cent at about 
10 8 Mq pc - 3 . If the first few Myr is the dominant binary-burning 
phase in the life-time of a cluster, current binary frequencies in 
clusters put constraints on the initial densities of open and globu- 
lar clusters, the respective ranges being indicated (Sec.[5]l. Lower 
panels: Deviations of best-fitting and model binary-fractions from 
the iV-body results. 



From the extrapolation it appears that the decrease of 
the binary-fraction flattens at higher densities. This may be 
a result of the choice of the functions describing the best- 
fitting parameters in Fig. (7) The bottle-neck here would be 
the — exp (— log 10 p ec i) behaviour chosen to describe the 5- 
parameter, that causes the slope of Q^J^'^ (S,t) to become 
constant beyond the maximum density reached in the inte- 
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grations. However, the flattening is already visible in the N- 
body based data between 10 3 and about a fewx 1O 5 M pc -3 
(solid lines in Fig. 110(1 . A physical interpretation is a very 
efficient energy generation due to the energy released in the 
more frequent interactions with hard binaries in dense sys- 
tems. This leads to strong cluster expansion in the densest 
configurations and therefore causes the binary evolution to 
end earlier than would be expected from extrapolating the 
binary fraction from low and intermediate densities to higher 
densities (Fig. H Sec. [3X2). 

As dissolution of binaries, and therefore the reduction of 
the binary-fraction, slows down with time (see the shrinking 
separation between the three lines in Fig. I10p . the evolu- 
tion of the binary population might halt already some Myr 
after the beginning of stimulated evolution, i.e. the binary 
pop ulation becomes frozen-in due to cluster expansion (see 
also Kroupa 1995a; Duchene et al.lll999l ; lFregeau et al . 2009; 
iParker et al.l 120091 ). This should be the case if roughly all 
soft binaries have dissolved. Especially for densities beyond 
« 10 s Mq pc -3 the binary fraction remains roughly con- 
stant already after 1 Myr, suggesting that binary-burning 
is so rapid that it has already reached an equilibrium sit- 
uation before 1 Myr is completed. If this first few Myr of 
evolution is the only or at least dominant part in the life of 
a cluster to alter the binary population, and if all clusters 
indeed start with /(, = 1, present-day observed binary popu- 
lations put constraints on the initial condition s under which 
the considered clusters should have fo rmed (|Kroupa et al.l 
ll999l ; lKroupa1l2000l ; |Parker et al.ll2009l ). 

Open clusters have global present-day bi nary- fractions 
in th e intermediate range (30 — 70 per cent, ISollima et al.l 
|2010| ) and should thus have formed with densities of 10 2 < 
Pod < 10 5 A / / Q pc -3 (see Fig. [TDJ. The average open clus - 
ter according to the catalogue by iPiskunov et alj (|2007l ) 
has r h fa 5.5 pcQ and M cc i = 2 x 10 3 M Q within its tidal 
boundary. The present-day mean density is thus about 
1.5Mq pc -3 only. None of the models end with such a 
low density at the end of the integration and our models 
can thus only be extrapolated. However, clusters in general 
can expand considerably as a result of the dropping poten- 
tial when the residual-gas fr om star formation is expelled 
l|Baumgardt fc Kroupal 120071 ). Binary depletion may thus 
have occured in the significantly denser embedded phase of 
cluster evolution. 

Pres ent-day total binary -fractions of 10 — 50 per cent 
in GCs (jSollima et al.l [2007) suggest their formation in 
somewhat denser environments, fewxlO 3 — 10 s Mq pc -3 . 
T hese initial densit i es are compatible with those derived 
in iMarks fc Kroupal |2010l ). who find p e d ~ 2 x 10 4 - 3 x 
1O 7 M0 pc -3 for the birth densities of individual GCs. An 
average Galactic GC today has a half-mass r adius of ~ 4 pc 
and a total mass of 4 x 1O 4 M (|Harridll99d 2003 revision, 
assuming a mass-to-light ratio of 1.5), resulting in a stellar 
mass-density within rh of w 7OM0 pc -3 . Indeed some of the 
JV-body model clusters from the stated initial density-range 
for GCs reach densities comparable to this value after 5 Myr 
of evolution (Fig. [4] lower panels). The density in a few ini- 
tially compact clusters that start with rh = 0.1 and 0.3 pc 



6 Using a conversion factor 1.85 valid for a Plummer model to 
calculate r\y from the r c -value that is given in the catalogue. 



decreases to « 1O 2 M0 pc -3 and below at the end of the in- 
tegration. The final densities in some rh = 0.8 pc models are 
also in the interesting density-range, but the binary-fraction 
remains too high in order to be a probable progenitor of a 
GC. 

The results can however also only be extrapolated 
to GCs, since none of our models starts or ends with a 
mass comparable to a GC. Additionaly G Cs will have been 
even more massive w hen they were born (|Marks k, KroupH 
l20ld ; IConrovl I 201 ll ). thus containing more massive stars 
and effects of stellar evolution might become important, too 
|lvanova et al.| [2005). It is also possible that during the life- 
time of a cluster further periods of binary burning occur. 
In particular the long-lived GCs can again reach densities 
sufficiently high to dissolve more binaries as they try to go 
into core-collapse. This would relax the need for very high 
initial densities. This has to be investigated in iV-body ex- 
periments that follow the dynamical evolution of a binary 
population for a longer time in higher-mass clusters. This 
is, however, not readily possible with present-day soft- and 
hardware. 



6 SUMMARY & OUTLOOK 

JV-body computations of the evolution over the first 5 Myr of 
binary populations in star clusters with initially 100 per cent 
binaries (similar to Oh et al., in prep.) are performed. Ini- 
tial orbit al parameters fo r the binaries are randomly selected 
from the lKroupal (|1995bf ) proto-binary birth period distribu- 
tion, masses for the components of binaries with primary- 
star masses < 5 Afp are selecte d randomly from the canon- 
ical stellar IMF (|Kroupall200lh . more massive binaries are 
paired to have a high mass-ratio, and eccentricities are se- 
lected from a thermal distribution. Before the integrations 
are started pre-main sequence eiqenevoluti on has been ap - 
plied to the so selected orbital parameters l| Kroupal [T995b) 
to obtain a realistic initial binary population. 

The integrations confirm earlier results that stimulated 
evolution first dissolves wide (low-energy) binaries due to 
their large cross-sectio n for intera ctions with other cluster 
members. We confirm iHutl |l983)'s finding that for fixed 
energy, binaries are more likely to be dissolved if they 
have a large reduced-mass, fi, but we also find that this 
p-dependence has a negligible effect on the results. 

The time-scale over which a population of primordial bi- 
naries in the present models (Tab. [I) evolves is found to be 
the crossing-time, i.e. the evolution of the whole binary pop- 
ulation in this set of computations can be well described by 
the initial cluster density, p cc i (t C r oc p~^' 5 )- Thus, initially 
denser clusters achieve lower binary-fractions than initially 
more extended configurations of the same age (the density- 
age degeneracy). How well this finding can be extrapolated 
to higher cluster masses needs to be researched in the future. 

Models which have different crossing- and relaxation- 
times but the same initial velocity dispersion, which sets the 
location of the hard-soft boundary for binaries, are found 
to have similar binary-fractions after a given number of 
relaxation-times. This is at present not fully understood, 
but is maybe related to the velocity dispersion in a clus- 
ter evolving on the energy-equipartition time-scale. If true 
then this suggests a close coupling of the binary population 
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and its host cluster. However, the early evolution of the ve- 
locity dispersion may rather be driven by stellar evolution 
and binary-burning than relaxation effects. The long-term 
two-body relaxation driven cluster evolution is not studied 
here, but the global properties of the binary population are 
mostly frozen by 5 Myr due to cluster expansion. 

Since the evolution of the same initial binary population 
in a cluster with a given density is unique, a stellar dynam- 
ical op erator, fijy" 1 ^ 11 (~t), can be introduced l|Kroupall2002l . 
l2008al ) which acts on the initial orbital-parameter binary 
distribution function, Dm, such that the resulting distribu- 
tion, T> Bcl ' rh (t), in the iV-body integrations is obtained. 
The operator f2 ( j„^ ;1 ' rh (t) is here quantified for the first time 
in terms of the initial cluster density and time for stimulated 
evolution. Therewith an analytical tool to efficiently calcu- 
late the evolved binary properties of a star cluster, given 
the above properties, is obtained. This recipe allows the ex- 
traction of the resulting energy-, period-, mass ratio and 
eccentricity distributions, ready for comparison with obser- 
vations. 

Assuming the initial binary population to have invariant 
properties and that the first few Myr of cluster evolution is 
dominant in changing a clusters' binary population the cur- 
rently observed binary properties of long-lived star clusters 
put constraints on the density at the time of their forma- 
tion. Open clusters which have observed binary frequencies 
of about 30 — 70 per cent should, according to these re- 
sults, have formed with densities fa 10 2 — 2 x 10 5 Mq pc -3 . 
Galactic globular clusters have binary- fractions of ~ 10 — 50 
per cent and must thus have formed somewhat denser, 
fewxlO 3 < peel < 10 8 Mq pc -3 , in excellent agreement wit h 
independently found constraints ( M arks fe Kroupal |2010| ). 
Although this work shows that cluster density is important 
in determining the binary-fraction in low-mass clusters, the 
constraints for GCs should be taken cautiously since none 
of the models have a mass comparable to present-day GCs. 
Thus the results can only be extrapolated. Additionaly the 
long-lived GCs may eventually undergo further periods of 
binary-burning through core-collapse episodes during the 
course of a Hubble-time and will have lost some systems 
due to evaporation and ejections. 

Applying the here developed tool to observed binary- 
populations in individual young star clusters will result in 
meaningful constraints for the mass and size of the star for- 
mation event from which the star cluster originated. The re- 
sults w ill be true only if the framework described in lKroupal 
|l995al lbllch is a valid description of physical reality, i.e. if all 
systems are born with the same BDFs, ft, = 1 and subse- 
quently undergo eigenevolution. However, the evolution of a 
binary population with the same relative occupancy of or- 
bits (i.e. the same form of t he period BDF) is inde pendent 
of the initial binary-fraction l|Kaczmarek et al1l201ll ) so that 
the present results and methods are expected to be applica- 
ble even if the initial binary fraction is less than unity. 

Additionaly, since all stars and binaries originate from 
discrete star formation events, galaxy-wide stellar popula- 
tions are the sum over all stars and binaries from all such 
events in their respectiv e host galaxy. In a follow-up paper 
|Marks fe Kroupap201ll ') we use this notion to calculate 
composite stellar populations of galactic fields and predict 
their properties (Dynamical Population Synthesis). This is 
an approach similar to calculating the integrated galactic 



IMF of stars in whole galaxies l|Kroupa fc Weidner! liool ; 
IWeidner fc Kroupall2005l ). 
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